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Detecting the limit cycles for a class of Hamiltonian systems 
under thirteen-order perturbed terms 

Gheorghe Tigan * 


Abstract 

It this paper we study a class of perturbed Hamiltonian systems under perturbations of thirteen 
order in order to detect the number of limit cycles which bifurcate from some periodic orbits of 
the unperturbed Hamiltonian system. The system has been previously studied in m, cni, m 
We observe in the present work that the system under perturbations of thirteen order can have 
more limit cycles than under perturbations of five cni, respectively nine im order but we have 
not identified more limit cycles than under perturbations of seven order. 
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1 Introduction 


Consider the following perturbed Hamiltonian system: 

f X = Ay {ahx"^ — by^ + l) + £x{ux" + vy^ — bj^x^y^ — ux"^ — A), , . 

\ y = Ax — aby"^ “ l) + £y{ux^ + uy" + bx^y^ — vy"^ — A) 

where /r + /3 = n, 0 < a < b < 1,0 < e<IC 1, it, w, A are the real parameters and n = 2k,k 
an integer positive. Studying the existence, number and distribution of limit cycles in a system of 
polynomial differential equations is an open problem even for a 2-degree polynomial and it is known 
as the Hilbert's 16th problem. 

The system of differential equations: 

J X = y {l + x^ — ay"^) + ex{mx^ + ny^ — A) 

\ y = —X (l — cx"^ + y^) + £y{mx‘^ + ny^ — A) ^ 

has been discussed in Q. E and the system 


J X = y (l — cy^) + £x{mx'^ + ny^ — A) 

^ 2 / = —X (l — ax^) -I- £y{mx‘^ + ny"^ — A) ^ 

in [7]. It has been shown that each of the two systems can have at least 11 limit cycles. 

The higher order perturbations of the system @, have been recently studied in 12], m 
The following result is reported in |2]. 
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Theorem 1.1 Consider the perturbed Hamiltonian system 


X =-—+ P{x,y,a),y = —+Q{x,y,a) (4) 

Assume that P{x,y,Q) = Q{x,y,Q) = 0, and the unperturbed system exhibits a center. Denote by 
the closed curves H{x,y) = h surrounding this center, and by T^{D) the subset bounded by and 
containing the center. Suppose that as h increases the diameter of the set r^{D) increases. If there 
exists ho such that function 


A{h)= J [P^^{x,y,0)+ Qy^{x,y,0)]dxdy (5) 

r'*(D) 

satisfies A[ho) = t),A'{ho) 7^ 0,aA'{ho) < 0, respectively > 0, then the perturbed system 0 has only 
one stable, respectively unstable limit cycle nearby F^”, for a very small. //F^ shrinks as h increases, 
the stability of the limit cycle is opposite to the above cases. If A{h) 0, then the system w has no 
limit cycle. 

The integral A{h) is called the Abelian integral If the Hamiltonian system 0 is perturbed in 
the form: 


f Ht) =-^+sx{pix,y) - X), 

I i'(^) = 11 + y) - ^). 

where p(0,0) = 9(0,0) = 0, then, by the above Theorem ll.il from A{h) = 0, we get: 


A = \{h) = 


J f(x,y)dxdy 
r'*(D) 

2 f dxdy 
r'*(D) 


( 7 ) 


with f(x, y) = a;p(. + ypy +p + q. This function X{h) is called the detection function of the system 0. 
Using the detection function X{h) we get by Theorem 1 1.1 1 the following result : 


Proposition 1.1 a) If (/iq, A(/io)) is an intersecting point of a line A = Aq and the detection curve 
A = X{h), with A'(/io) > 0(< 0), then the system 0) has only one stable (unstable) limit cycle nearby 
F^“; b) If the line A = Aq and the detection curve A = X{h) do not intersect each other, then the 
system W has no limit cycle. If the r^(D) shrinks as h increases, the stability of the limit cycle is 
opposite to the above cases. 


This paper is organized as follows. In Section 2, we recall the results regarding the unperturbed 
system. In Section 3, using the analytical expressions of the detection functions and a Computer 
Algebra System we numerically compute the detection functions. Based on these data we can illustrate 
the distribution of the limit cycles. 


2 The behavior of the unperturbed system 

The unperturbed system ® of the system o is: 
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(8) 


J X = Ay (—by'^ + abx'^ + l) 

\ y = Ax (ax^ — aby'^ — l) 

The Hamilton function defining this system is: 

H (x, y) = —(ax^ + by'^) + 2abx'^y^ + 2 + y^) = h (9) 

The function H has nine finite singular points: 

0 ( 0 , 0 ), v'6(l-6a)(l + 6)) , A , ^b{l - ba) { 1 + b)) , 

^3 ^(T^, (1 - M (1 + b)) ,A4 , 

Hi(0, ^),S 2 ( 0 , -yT)>Oi(^,0) and C 2 {-^,,Q) ■ 

Computing the eigenvalues of the associated linear system at each singular point we get that the 
points O, Hi, ^ 2 , ^ 3 , A 4 are centers, while i?i, H 2 , Ci,C 2 are saddle points. 

The values of the Hamilton function H at the singular points are respectively: H {Ai) = i = 

1 — A, H {Bk) = {Ck) = = 1,2. Because 0<a<6<lwe get that: H {O) < H (Hi) < 

i?(Ci) <i7(Hi). 

In polar coordinates, x = r cos 0,y = rsin0, the system becomes: 

f = —r^p'{9), 0 = —1 + Ep{9) (10) 

and the Hamilton function gets: 

H{r,9) =—r'^plff) + 2 E = h, (11) 

where 

p{9) = a cos"^ 9 + &sin^ 9 — 2a6cos^ 0sin^ 9. (12) 


Remark 


1 


The equilibrium points Hi,H 2 ,H 3 ,H 4 lie on the lines d± : 9 = iarccos 


b-\-ba 

a-\-b-\-2ab 


Theorem 2.1 As h varies on the real line, the level curves H{x,y) = h can be classified, as follows 

m 

1 . Ti : —00 < h < 0, which corresponds to an orbit that surrounds all critical points (FigWh). 

2. r2 UTi : 0 < h < ^, corresponding to an orbit (T^) that surrounds only the point O and a curve 
of type (Tf) , (Fig^,a). 

3. Vlf ^ < h < ^, which corresponds to two symmetric orbits that do not cross the Oy axis, but 
encircle the rest of the critical points. If h = ^ we get four heteroclinic orbits connecting the 
critical points Bi and B 2 (Fig\^,a). 

r 4 : i < h < ’ corresponding to four orbits that surround respectively the points 

Hi,i = 1 — 4. If h = ^ we have four homoclinic orbits, namely, two homoclinic to Ci and two 
homoclinic to C 2 (Fig\^,a). Note that as h increases, the curves r5*,r3 and T\ shrink, while 
T 2 extends. 
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Figure 1: Orbit of type a) Li (left) b) L 2 and Li (right) 



Figure 2: a) Four heteroclinic orbits connecting two critical points Bi, B 2 (left) b) Two orbits of 
type L 3 (right) 


3 Numerical explorations 


In this section we numerically compute the detection curves and point out the distribution of the limit 
cycles. The four detection functions can be computed numerically, and for a given h, they depend on 
u and V, (see tables 1-4). On the other hand, for two given values of u and v, the detection curves 
can be plotted on the (h, A)-plane, as illustrated in FigQIEl By the Pronosition ll .1 l aud the detection 
function graphs, the existence, number and distribution of limit cycles can then be illustrated. We 
consider here the case n = 12. 

From HU, we get 


ru2 = r±(6», h) 


1±^1- hp{0) 

pie) 


(13) 


and from 9 = —\ + r^p{9) = 0, we have: 


9i{h) = I arccos 
e 2 {h) = ^arccos 


(^b — a + 2^a?b'^ — ab {ab + 2ab) h 
[b — a — 2\JaP'b'^ — ab -\- [a + b + 2ab) 


/ b 2(2^) 

/ (^Q, b 2q5) 


Using the perturbation terms 

p (x, y) = X (ux"^ + uy" — b^^x^y^ — ux'^'\ ,q{x,y) = y (itx” -I- uy” -I- bx'^y^ — vy"^) 
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Figure 3: a) Four homoclinic orbits connecting the critical points Ci and C 2 (left) b) Four orbits of 
type L 4 (right) 


we have ^ ^ ayd’e'^ = (2 + n) (ux" + vy"^) — 3{ux^ + vy'^) — 2A. 

Hence, the four detection functions, corresponding to the four closed curves F^, j = 1 — 4, for the 
above perturbations are as follows: 


\,{h) 


f [(n + 2) (ux^ + vy’^) — 3{ux‘^ + vy"^)] dxdy 
r'jiD) 

2 f dxdy 
r^iD) 


(14) 


For a = 1/3,5 = 1/2, and n = 12, we get by itTHl the next four detection functions, in polar 
coordinates: 


2-17 

f (rl {0, h) g (9) - |r? {9, h) gi {9)) d9 
Ai(/i) = - — -, —00 < h <2, 

J {^1 h) dO 

0 


/ {rl {0, h) g (6») - f {9, h) gi (9)) d9 

Mh) = - --, 0 < h < 2, 

f r2 (9, h) d9 


(15) 


(16) 


62{h) 

f [{rl {d, h) - rl {9, h)) 5 ( 0 ) - | {rf {9, h) - rj {9, h)) gi ( 6 »)] d9 

X,{h) = - — -, 2<h<3. (17) 

/ (ri{9,h) -r2{9,h))d9 

-e2(h) 
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3</i<8.2, (18) 


Hh) = 


e2{h) 

f [(’’i (6', h) - rl {9, h)) g (6») - | (rf {9, h) - rl (6», h)) gi (6»)] d9 

ei(h) 

MM ’ 

/ {ri{9,h) -r2{9,h))d9 

9lih) 


where g (9) = u cos^^ 9 + v sin^^ 9, gi (9) = u cos^ 9 + v sin^ 9 and ri _2 {9, h) = (0, h). 

Using the expressions (P )l - (P|l and a Computer Algebra System we find the values of the detection 
functions Xi{h), i = 1 — 4, recorded in the following tables (1-4). Denote by p = lO^u and w = lO^w 


Table 1 

Values of the detection function Xi{h), for a = 1/3, b = 1/2, n = 12. 


h 

Xi{h) 

h 

Xi{h) 

h 

Xi{h) 

-2 

4.933p-f 1.373u; 

-1.9 

4.862p-t 1.352a; 

-1.8 

A.792p+ 1.332a; 

-1.7 

4.722p-t 1.311W 

-1.6 

4.653p-f 1.291a; 

-1.5 

4.584p+ 1.271a; 

-1.4 

4.516p-|- 1.251a; 

-1.3 

4.448p-f 1.231a; 

-1.2 

4.381p-b 1.212a; 

-1.1 

4.315p-f 1.192a; 

-1. 

4.249p-f 1.173a; 

-0.9 

4.183p-b 1.154a; 

-0.8 

4.118p-f 1.135a; 

-0.7 

4.054p-t 1.116a; 

-0.6 

3.990p-b 1.098a; 

-0.5 

3.926p-f 1.079a; 

-0.4 

3.863p-f 1.061a; 

-0.3 

3.801p-b 1.043a; 

-0.2 

3.739p-f 1.025a; 

-0.1 

3.678p-f 1.008a; 

0 

3.618p-f 0.9906a; 

0.1 

3.558p-b 0.9733a; 

0.2 

3.498p-f 0.9562a; 

0.3 

3.439p-f 0.9392a; 

0.4 

3.381p-b 0.9224a; 

0.5 

3.323p-f 0.9058a; 

0.6 

3.266p-f 0.8895a; 

0.7 

3.210p-b 0.8733a; 

0.8 

3.154p-t 0.8573a; 

0.9 

3.099p-f 0.8415a; 

1 . 

3MAp+ 0.8259a; 

1.1 

2.990p-f 0.8105a; 

1.2 

2.937p-f 0.7954a; 

1.3 

2.885p-b 0.7805a; 

1.4 

2.834p-t 0.7658a; 

1.5 

2.783p-f 0.7514a; 

1.6 

2.734p-h 0.7373a; 

1.7 

2.685p-f 0.7235a; 

1.8 

2.638p-f 0.7102a; 

1.9 

2.593p-b 0.6974a; 

2 . 

2.553p-f 0.6859a; 



Table 2 

Values of the detection function X 2 {h), for a = 1/3, b = 1/2, n = 12. 


h 

X2{h) 

h 

X2{h) 

0.01 

-0.001875 u-0.001876 v 

0.11 

-0.020732 u-0.02083 v 

0.21 

-0.039766 u-0.0401393 v 

0.31 

-0.058969 u-0.059813 v 

0.41 

-0.078305 u-0.079837 v 

0.51 

-0.097701 u-0.10015 v 

0.61 

-0.117021 u-0.120612 v 

0.71 

-0.136046 U-0.140967V 

0.81 

-0.154444 u-0.160774 v 

0.91 

-0.171734 u-0.179321 v 

1.01 

-0.187236 u-0.195488 v 

1.11 

-0.200019 u-0.207541 v 

1.21 

-0.208827 u-0.212808 v 

1.31 

-0.211988 u-0.207141 v 

1.41 

-0.207302 u-0.184002 v 

1.51 

-0.19189 u-0.132774 V 

1.61 

-0.162016 u-0.0354102 v 

1.71 

-0.11285 u-1-0.141107 V 

1.81 

-0.038213 u-l-0.465149 v 

1.91 

0.069632 u-f 1.1178 v 
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Table 3 

Values of the detection function X 3 {h), for a = 1/3,6= 1/2, n = 12. 


h 

Xsih) 

h 

X3{h) 

h 

X3{h) 

2 

3.0408p -f 0.8167W 

2.02 

3.0446p 0.8175a; 

2.04 

3.0459p -h 0.8177a; 

2.06 

3.0463p -f 0.8177W 

2.08 

3.0462p 0.8175a; 

2.1 

3.0457/O -h 0.8172a; 

2.12 

3.0449p -l-0.8168a; 

2.14 

3.0438p -h 0.8164a; 

2.16 

3.0425p -h 0.8159a; 

2.18 

3.0411p -h 0.8154a; 

2.2 

3.0396p 0.8148a; 

2.22 

3.0379p -h 0.8142a; 

2.24 

3.0361p -f 0.8136a; 

2.26 

3.0342p 0.8130a; 

2.28 

3.0323p -h 0.8123a; 

2.3 

3.0303p -f 0.8116a; 

2.32 

3.0283p 0.8110a; 

2.34 

3.0262p -h 0.8103a; 

2.36 

3.0240p -f 0.8096a; 

2.38 

3.0219p -h 0.8089a; 

2.4 

3.0197p -h 0.8082a; 

2.42 

3.0176p -f 0.8075a; 

2.44 

3.0154p -h 0.8068a; 

2.46 

3.0132p -h 0.8061a; 

2.48 

3.0110p -f 0.8054a; 

2.5 

3.0089p 0.8047a; 

2.52 

3.0068p -h 0.8040a; 

2.54 

3.0047p -f 0.8033a; 

2.56 

3.0026p 0.8027a; 

2.58 

3.0006p -h 0.8020a; 

2.6 

2.9986p -f 0.8014a; 

2.62 

2.9967p -h 0.8008a; 

2.64 

2.9949p -h 0.8002a; 

2.66 

2.9931p -f 0.7996a; 

2.68 

2.9915p -h 0.7991a; 

2.7 

2.9899p -h 0.7986a; 

2.72 

2.9884p -h 0.7981a; 

2.74 

2.9871p 0.7976a; 

2.76 

2.9859p 0.7972a; 

2.78 

2.9848p -f 0.7968a; 

2.8 

2.9839p 0.7965a; 

2.82 

2.9833p -h 0.7963a; 

2.84 

2.9828/0 -f 0.7961a; 

2.86 

2.9826p 0.7959a; 

2.88 

2.9827/0 -h 0.7959a; 

2.9 

2.9831p -f 0.7959a; 

2.92 

2.9840p 0.7961a; 

2.94 

2.9854/0 -h 0.7964a; 

2.96 

2.9876p -f 0.7970a; 

2.98 

2.9909p 0.7978a; 

3. 

2.9973p -h 0.7995a; 


Table 4 


Values of the detection function X 4 {h), for a = 1/3, b = 1/2, n = 12. 


h 

Xiih) 

h 

X4{h) 

h 

X4{h) 

3 

2.9973p-H0.7995a; 

3.04 

3.0074p-h0.8022a; 

3.08 

3.0118p-f0.8032a; 

3.12 

3.0140p-H0.8037a; 

3.16 

3.0149p-f0.8037a; 

3.2 

3.0147p-f0.8035a; 

3.24 

3.0138p-H0.8030a; 

3.28 

3.0121p-h0.8024a; 

3.32 

3.0098p-f0.8016a; 

3.36 

3.0070p-H0.8006a; 

3.4 

3.0038p-f0.7995a; 

3.44 

3.0001p+0.7982a; 

3.48 

2.9960p-l-0.7969a; 

3.52 

2.9915p-f0.7954a; 

3.56 

2.9867p-f0.7939a; 

3.6 

2.9816p-l-0.7922a; 

3.64 

2.9762p-h0.7905a; 

3.68 

2.9706p-f0.7887a; 

3.72 

2.9646p-H0.7868a; 

3.76 

2.9584p-h0.7849a; 

3.8 

2.9520/o-h0.7829a; 

3.84 

2.9453p-H0.7808a; 

3.88 

2.9384p-h0.7787a; 

3.92 

2.9313p-f0.7765a; 

3.96 

2.9240p-l-0.7742a; 

4. 

2.9165p-f0.7719a; 

4.2 

2.8764p-h0.7596a; 

4.4 

2.8322p-l-0.7461a; 

4.6 

2.7847p-f0.7317a; 

4.8 

2.7341p-h0.7165a; 

5. 

2.6808p-H0.7005a; 

5.2 

2.6250p-f0.6838a; 

5.4 

2.5671p+0.6665a; 

5.6 

2.5072p-H0.6487a; 

5.8 

2.4455p-f0.6304a; 

6. 

2.3820p-f0.6116a; 

6.2 

2.3171p-H0.5924a; 

6.4 

2.2507p-f0.5729a; 

6.6 

2.1830p-f0.5530a; 

6.8 

2.1141p-H0.5327a; 

7. 

2.0440p-f0.5122a; 

7.2 

1.9729p-f0.4914a; 

7.4 

1.9008p-H0.4704a; 

7.6 

1.8278p-f0.4491a; 

7.8 

1.7540p-f0.4277a; 

8. 

1.6794p-H0.4060a; 

8.2 

1.6041p-h0.3842a; 



From the tables (1-4) we have the four detection functions, three of which are shown in FigQlEl 
The values of A 2 are not plotted because they are too small in comparison with the other values. From 
ProDOsition ll.il and FigQIEl one gets: 
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Figure 4: Detection curve A4 of the system , for a = = ^,n = 12, u = 0.007 and v = —0.028. 



Figure 5: Detection curves A3 (left), Ai (right) of the system (^1, for a = ^,b = ^,n = 12,u = 0.007 
and V = —0.028 


Theorem 3.1 For a = ^,b = ^,n = 12,u = 0.007 and v = —0.028 and 0 < £ <C 1, we have the 
following distribution of the limit cycles: 

a) If 176.22 < A < 242.6, the system m has at least four limit cycles, one in the neighborhood of each 
orbit of type T\, Fig\^), 

b) If 242.6 < A < 286.76, the system 0) has at least five limit cycles, one in the neighborhood of 
each orbit of type T\ and Ff, Fig\^), 

c) If 286.76 < A < 288.49, the system m has at least nine limit cycles, two in the neighborhood 
of each orbit of type Fg, and one in the neighborhood of each orbit of type T\, Fg, Fig\^), 

d) If 288.49 < A < 288.92, the system Q) has at least eleven limit cycles, two in the neighbor¬ 
hood of each orbit of type T\ and one in the neighborhood of each orbit of type 

e) If 288.92 < A < 289.99, the system m has at least thirteen limit cycles, two in the neighbor¬ 
hood of each orbit of type T\, Fg and one in the neighborhood of the orbit of type Ff, Fig^), 

f) If 289.99 < A < 290.82, the system ^ has at least nine limit cycles, two in the neighborhood 

of each orbit of type F^ and one in the neighborhood of the orbit of type r^, Fig^) . 

From FigEEl we could obtain more results but we listed only the important cases. 
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Figure 6: Distribution diagram corresponding to: a) four (left), b) five (middle), c) nine (right) limit 
cycles of the system o 



Figure 7: Distribution diagram corresponding to: a) eleven (left), b) thirteen (middle), c) nine (right) 
limit cycles of the system o 


4 Conclusion 

The system Q corresponding to a = ^,b = = 12, u = 0.007, u = —0.028,0 < e <C 1, and 

288.49 < A < 288.92, has at least eleven limit cycles while for 288.92 < A < 289.99 it has at least 
thirteen limit cycles. The Abelian integral method was used. Through numerical explorations we have 
drawn the shape of the graphs of the detection functions, from which we determined the number and 
distribution of limit cycles. A natural question is: Is it possible to obtain more limit cycles for higher 
perturbations? 
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